# Paper: Is it still the economy? Economic voting in polarized politics
# Individual-level simulation
# Replication: Figures S10 and S11 (Supplementary Material)
remove(list = ls())

library(conflicted)
library(tidyverse)
library(dplyr)
conflict_prefer("filter","dplyr")
library(ggplot2)
library(here) # When we have an R-project, it reads files without setting an directory
conflict_prefer("here", "here")
library(nnet)
library(ggeffects)
library(MASS)
conflict_prefer("select","dplyr")

load(here("Final_Paper", "R&R", "CountyLevel.RData"))

polar_data <- model_data %>% 
  select(year, elitepolar_house) %>% 
  group_by(year) %>% 
  summarise(pol_median_dist1 = mean(elitepolar_house, na.rm = T)) %>% 
  ungroup()

load(here("Final_Paper", "R&R", "IndLevel.RData"))

model_data <- left_join(model_data, polar_data, by = "year")
glimpse(model_data)


summary(ml_pol <- lm(percep_polar ~ pol_median_dist1*ptd_id3 + 
                       female + black + hispanic + asian + 
                       other_race + log_age + college + income + state_fips, 
                     model_data))

model_pol <- ggeffect(ml_pol, terms = c("pol_median_dist1[0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95]", 
                                        "ptd_id3"), 
                      ci.lvl = 0.95)
glimpse(model_pol)
table(model_pol$group)

f2 <- ggplot(model_pol, aes(x = x, y = predicted, group = group)) +
  geom_line(aes(linetype = group, color = group), size = .75) + 
  theme_light() +
  scale_linetype_manual(name = "Partisanship:",
                        values = c("solid", "longdash", "dashed")) +
  scale_color_manual(name = "Partisanship:",
                     values = c("darkblue", "black", "darkred")) +
  geom_ribbon(aes(ymin = conf.low, ymax = conf.high, fill = group), alpha = .1) +
  scale_fill_manual(name = "Partisanship:",
                    values = c("blue", "gray30","red")) +  # #7fb9b4 and #bc7c86
  ggtitle(NULL) +
  xlab("Elite Polarization in Congress") +
  ylab("Perceived Ideological Gap") +
  #ggtitle("(b) Relationship between Observed and Perceived Polarization") +
  #labs(linetype = "Party Id.") +
  theme(legend.position="bottom",
        axis.title = element_text(size = 13),
        axis.text = element_text(size = 12),
        legend.title = element_text(size = 13),
        title = element_text(size = 10),
        legend.text = element_text(size = 12))
f2

summary(ml_pol2 <- lm(percep_polar ~ pol_median_dist1*ptd_id5 + 
                        female + black + hispanic + asian + 
                        other_race + log_age + college + income + state_fips, 
                      model_data))

model_pol2 <- ggeffect(ml_pol2, terms = c("pol_median_dist1[0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95]", 
                                          "ptd_id5"), 
                       ci.lvl = 0.95)
glimpse(model_pol2)
table(model_pol2$group)
f3 <- ggplot(model_pol2, aes(x = x, y = predicted, group = group)) +
  geom_line(aes(linetype = group, color = group), size = .75) + 
  theme_light() +
  scale_linetype_manual(name = "Partisanship:",
                        values = c("solid", "longdash", "dashed", 
                                   "twodash","dotdash")) +
  scale_color_manual(name = "Partisanship:",
                     values = c("darkblue", "blue", "black",
                                "red", "darkred")) +
  geom_ribbon(aes(ymin = conf.low, ymax = conf.high, fill = group), alpha = .1) +
  scale_fill_manual(name = "Partisanship:",
                    values = c("blue", "#7fa1b9", "gray30",
                               "#ea4e67", "red")) +  # #7fb9b4 and #bc7c86
  ggtitle(NULL) +
  xlab("Elite Polarization") +
  ylab("Perceived Ideological Distance") +
  #labs(linetype = "Party Id.") +
  theme(legend.position="bottom",
        axis.title = element_text(size = 13),
        axis.text = element_text(size = 11),
        legend.title = element_text(size = 13),
        legend.text = element_text(size = 11))
f3


observed_data <- polar_data %>% 
  group_by(year) %>% 
  #mutate(pol_median_dist1 = pol_median_dist1/2) %>% 
  summarise(gap = mean(pol_median_dist1, na.rm = T)) %>% 
  ungroup() %>% 
  mutate(ptd = "Elite Polar")
summary(observed_data$gap)

rep_data <- model_data %>% filter(ptd_id3 == "Rep")
annualrep_data <- rep_data %>% 
  group_by(year) %>% 
  summarise(gap = mean(percep_polar, na.rm = T)) %>% 
  ungroup() %>% 
  mutate(ptd = "Rep")

dem_data <- model_data %>% filter(ptd_id3 == "Dem")
annualdem_data <- dem_data %>% 
  group_by(year) %>% 
  summarise(gap = mean(percep_polar, na.rm = T)) %>% 
  ungroup() %>% 
  mutate(ptd = "Dem")

ind_data <- model_data %>% filter(ptd_id3 == "Ind")
annualind_data <- ind_data %>% 
  group_by(year) %>% 
  summarise(gap = mean(percep_polar, na.rm = T)) %>% 
  ungroup() %>% 
  mutate(ptd = "Ind")

annual_data <- model_data %>% 
  group_by(year) %>% 
  summarise(gap = mean(percep_polar, na.rm = T)) %>% 
  ungroup() %>% 
  mutate(ptd = "All Voters")
glimpse(annual_data)

gaps_data <- rbind(observed_data, annual_data, annualrep_data, 
                   annualdem_data, annualind_data)
glimpse(gaps_data)

glimpse(annualrep_data)
ts_data <- left_join(annualrep_data, annualdem_data, by = "year")
glimpse(annualind_data)
ts_data <- left_join(ts_data, annualind_data, by = "year")

glimpse(gaps_data)
gaps_data <- gaps_data %>% 
  mutate(ptd = factor(ptd,
                      levels = c("Elite Polar",
                                 "All Voters",
                                 "Dem",
                                 "Ind", 
                                 "Rep"))) 
table(gaps_data$ptd)


gaps_data <- gaps_data %>% 
  mutate(gap = ifelse(ptd == "Elite Polar", gap/1.5, gap))


gaps <-  gaps_data %>% 
  ggplot(aes(x = year, y = gap, color = ptd, shape = ptd, fill = ptd, 
             linetype = ptd)) +
  geom_point(size = 2.7) + 
  geom_line(size = 1.2) + 
  theme_light() +
  scale_color_manual(name = NULL,
                     values = c("black", "orange", "blue", "gray45", "red"),
                     labels = c("Elite", "All Voters", "Dem", "Ind", "Rep")) +
  scale_linetype_manual(name = NULL,
                        values = c("solid", "longdash", "dashed", "twodash", "dotdash"),
                        labels = c("Elite", "All Voters", "Dem", "Ind", "Rep")) +
  scale_shape_manual(name = NULL,
                     values = c(21, 22, 24, 23, 25),
                     labels = c("Elite", "All Voters", "Dem", "Ind", "Rep")) +
  scale_fill_manual(name = NULL,
                    values = c("black", "orange", "blue", "gray45", "red"),
                    labels = c("Elite", "All Voters", "Dem", "Ind", "Rep")) +
  ylab("Perceived Ideological Gap") + xlab(NULL) +
  scale_x_continuous(breaks = c(1992, 1996, 2000, 2004, 2008, 2012, 2016, 2020)) +
  theme(axis.text.x = element_text(hjust = 0.5, size = 11),
        axis.text.y = element_text(size = 11),
        axis.title.y = element_text(size = 14),
        axis.title.x = element_text(size = 14),
        title = element_text(size = 10),
        legend.position = "bottom",
        legend.direction = "horizontal",
        legend.text = element_text(size = 14))
gaps

gaps <- gaps + scale_y_continuous(sec.axis = sec_axis(~.*1.5, 
                                                      name = "Observed Elite Polarization"))

gaps
